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Abstract. We study the dynamics of long-wavelength fluctuations in one- 
dimensional (ID) many-particle systems as described by self-consistent mode- 
coupling theory. The corresponding nonlinear integro-differential equations 
for the relevant correlators are solved analytically and checked numerically. 
In particular, we find that the memory functions exhibit a power-law decay 
accompanied by relatively fast oscillations. Furthermore, the scaling behaviour 
and, correspondingly, the universality class depends on the order of the leading 
nonlinear term. In the cubic case, both viscosity and thermal conductivity diverge 
in the thermodynamic limit. In the quartic faster decay of the memory 

functions leads to a finite viscosity, while thermal conductivity exhibits an even 
faster divergence. Finally, our analysis puts on a more firm basis the previously 
conjectured connection between anomalous heat conductivity and anomalous 
diffusion. 
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1. Introduction 

Numerical and theoretical studies of transport in reduced spatial dimensions (1 and 
2D) have revealed that unexpected properties of matter may emerge with respect 
to the standard features of 3D systems. This fact, combined with the possibility 
of designing and manipulating matter at the nanoscale, unveals the possibility of 
exploring new physical phenomena and applications. Preliminary indications about 
anomalous hydrodynamic behaviors were discovered in the context of non-equilibrium 
statistical mechanics by the appearance of long-time tails in fluid models [TJ G] . In 
fact, there are many instances showing that microscopic models of hydrodynamics 
can be affected by serious pathologies. For the sake of an example, we mention the 
anomalous scaling of Rayleigh and Brillouin peak widths in the hydrodynamic limit 
detected in the ID Lennard- Jones fluid [21 0]. Another instance is the divergence 
of viscosity in cellular automata fluids [5] that may however be due their peculiar 
features. 
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On the one hand, the breakdown of phenomenological constitutive laws of 
hydrodynamics, yielding ill-defined transport coefficients, could be viewed as a 
drawback for the adopted models and theoretical approaches. On the other hand, it is 
nowadays widely recognized that anomalous transport is a documented effect rather 
than an accident due to the oversimplification of theoretical models. Typical examples 
are single-filing systems, where particle diffusion does not follow Fick's law [6], and 
the enhancement of vibrational energy transmission in quasi-lZ) systems, such as 
polymers [7] or individual carbon nanotubes [8]. 

Within this general context, one of the issues that attracted a renewed interest in 
the last decade is the problem of anomalous heat conduction in low-dimensional many- 
particle systems [3J[T0]. In this case, the anomalous features amount to the divergence 
of the finite-size heat conductivity k(L) oc L a in the limit L — > oo and, correspondingly, 
to a nonintegrable decay of equilibrium correlations of the energy current (the Green- 
Kubo integrand), (J(t)J(O)) oc t~^ a \0 < a < 1) for long times t — > oo. Numerical 
studies [TT] and theoretical arguments [T3] indicate that anomalies should occur 
generically in 1 and 2D, whenever momentum is conserved. In particular, a great 
deal of studies have been devoted to the prediction of the exponent a for a twofold 
reason. On a basic ground, one is interested in identifying the universality classes and 
the relevant symmetries possibly identifying them. From a practical point of view, an 
accurate prediction of the scaling behaviour allows determining heat conductivity in 
finite systems - a crucial issue for the comparison with experiments. 

Despite many efforts, the theoretical scenario is still controversial. The first 
attempt of tackling the problem was based on arguments derived from the Mode- 
Coupling Theory (MCT). This approach proved quite powerful for the estimate of 
long-time tails in fluids [I] and for describing the glass transition [T3]. In ID models 
[14] . the MCT prediction was a = 2/5 [TDl [15], very close to the numerical values 

This estimate was later criticized as inconsistent in [12], where renormalization 
group arguments were instead shown to yield a — 1/3. Nevertheless, the 2/5 value has 
been later derived both from a kinetic-theory calculation for the quartic Fermi-Pasta- 
Ulam (FPU-/?) model [TB] and from a solution of the MCT by means of an ad hoc 
Ansatz [T7]. It was thereby conjectured [T7] that 2/5 is found for a purely longitudinal 
dynamics, while a crossover towards 1/3 is to be observed only in the presence of a 
coupling to transversal motion. Unfortunately, the accuracy of numerical simulations 
is generally insufficient to disentangle the whole picture. Without pretending to give a 
full account of the various numerical results (see |18] for a recent account), it is worth 
mentioning that, while for some models like the so called hard-point gas [TH] the data 
is compatible with the prediction a — 1/3 [20] , some other cases show huge deviations. 
For example, in the purely quartic FPU model, a is definitely larger than 1/3 and 
definitely closer to 2/5 [21]. The situation is even more controversial in 2D, where 
MCT predicts a logarithmic divergence 22J that is not yet confirmed by numerical 
simulations. 

For all the above reasons, we decided to reconsider MCT to clarify the origin 
of the reported discrepancies. Typically, MCT is implemented by imposing an ad- 
hoc Ansatz for the dependence of the hydrodynamic propagator and of the memory 
kernel on the space and time variables. In this paper, following [23], we find an 
exact self-consistent solution of the MCT equations, without any a priori assumption 
on the form of the solution. The resulting scaling behavour differs from that one 
proposed in [17] . Moreover, the asymptotic behaviour turns out to depend on the 
order of the leading nonlinearity in the interaction potential. Indeed, we find that 
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cubic and quartic nonlinearities are characterized by a = 1/3, and 1/2, respectively. 
Remarkably, an even more striking difference is exhibited by viscosity which appears 
to diverge in the former case, while is finite in the latter one. 

The structure of the paper is as follows. In Section [5] we review the mode- 
coupling equations and discuss their physical meaning. Section [3] is devoted to the 
solution when the leading nonlinearity is cubic, while in Section [4] we discuss how the 
solutions change in the presence of a leading quartic nonlinearity. In Section [5j the 
knowledge of the mode-coupling solutions is exploited to determine the value of the 
scaling exponent a, by estimating the heat current correlation function. Final remarks 
and future perspectives of this approach are reported in the conclusions. 

2. MCT equations 

The introduction of suitable stochastic equations ruling the dynamics of the relevant 
variables is rather customary to describe the relaxation close to equilibrium [24| . The 
idea consists in concentrating on the effective motion of the "slow" observables, i.e. in 
changing level of description. The general strategy amounts to projecting the original 
equations onto a suitable subspace and eventually leads to linear non-Markovian 
equations. The memory term determines the relaxation properties. Whenever a sharp 
separation of time scales exist, the equations reduce to their Markovian limit. 

Let us consider the simplest one-dimensional version of the self-consistent MCT. 
For dcfmitcncss, we consider a Hamiltonian chain of oscillators, interacting through 
a generic nearest-neighbour coupling potential V, whose Taylor expansion around 
equilibrium is written as 

V(y) = \y 2 + lg 3 y 3 + + ... (1) 

The variable y is a short-hand notation for the difference of nearest neighbour 
displacement fields — u.;), with the integer index i labelling the lattice sites. 

Besides the conservation of total energy, such a kind of interaction implies also the 
conservation of total momentum due to space-translation invariance. 

The main observable we are interested in is the normalized correlator G{q, t) = 
(Q*(q, t)Q(q, 0))/(|Q(q)| 2 ), where Q(q, t) is the Fourier transform of the displacement 
field Ui(t). By assuming periodic boundary conditions for a chain made of N sites, 
the wavenumber is given by q = 2irk/N, with —N/2 + 1 < k < N/2. Notice also that 
G(q, t) = G(—q, t). We simplify the notation by setting to unity the particle mass, the 
lattice spacing and the bare sound velocity. The equations for the correlator G(q, t) 
then read [251 CDS 

G(q,t)+e [ T(q,t- s)G(q,s)ds + uj 2 (q)G(q,t) =0 (2) 
Jo 

where the memory kernel T(q,t) is proportional to (J-(q,t)J-(q,0)) with J-(q) being 
the nonlinear part of the force between particles. Equations ([2]) must be solved with 
the initial conditions G(q, 0) = 1 and G(q, 0) = 0. Equations (|2|) are exact and they 
are derived within the well-known Mori-Zwanzig projection approach [2"4] , 

The mode-coupling approach basically amounts to replacing the exact memory 
function T with an approximate one, where higher-orders correlators are written 
in terms of the G(q,t). This yield a closed system of nonlinear integro-differential 
equations. For potentials like ([T]) this has been worked out in detail in references 
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|251 I15j . Both the coupling constant e and the frequency w(q) are temperature- 
dependent input parameters, which should be computed independently by numerical 
simulations or approximate analytical estimates [251 I15j . For the aims of the 
present work, we may restrict ourselves to considering their bare values, obtained 
in the harmonic approximation. In the adopted dimensionless units they read 
e = ig^ksT /2ir and u)(q) = 2| sin||. Of course, the actual renormalized values are 
needed for a quantitative comparison with specific models. 

To understand the physical interpretation of equation (0), let us first assume that 
a Markovian approximation holds i.e. we can replace the memory function by a Dirac 
delta. In the small wavenumber limit q — > they reduce to 

G(q,t)+ r ,q 2 G(q,t)+c 2 q 2 G(q,t) = (3) 

where, for the sake of clarity, we have reintroduced the sound velocity c. Notice 
that the validity of ((3|) rests on the fact that T exhibits a fast decay and that its 
time integral is finite. This equation describes the response of an elastic string at 
finite temperature as predicted by linear elasticity theory |26j for the macroscopic 
displacement field u(x, t) 

u-c 2 dlu-r)d 2 u = . (4) 

Here r\ is a suitable viscosity coefficient and describes the internal irreversible processes 
of the elastic string (i.e. the absorption of "sound" waves). We anticipate that one of 
our results is precisely that (generically) the above approximation does not hold and 
that long-range memory effects must be retained |27j . 



3. Cubic nonlinearities 



In the generic case in which (73 in equation |T]) is different from zero, the lowest-order 
mode coupling approximation of the memory kernel turns out to be |25l 115] 

T(q,t)=u, 2 (q)- £ G(p,t)G(p',t) . (5) 

p+p' — q— 0,±7T 

Equation @ along with are now a closed system that has to be solved self- 
consistently. The integral terms contain products of the form GGG and we thus 
refer to this case as to the cubic nonlinearity. 



3.1. Analytical results 

Direct numerical simulations [15) indicate that for small wavenumbers nonlinear and 
nonlocal losses in equation (|2|) are small compared with the oscillatory terms. This 
suggests splitting the G dynamics into phase and amplitude evolution, 

G(q,t) = C(q,0e M(<z) * +c.c. (6) 

As already shown in [23], this Ansatz allows reducing the dynamics to a first order 
equation for C(q,t) and, in turn, determining the scaling behavior of the correlators. 
Let us start by rearranging the expression for r(g, t). By inserting equation ^ into 
(O, we obtain 



F(q,t) = uj 2 {q) / dp C{p,t)C{q-p,t)e l ^ + ^ q -^ t (7) 

J — 7T 

+ C{p,t)C*(q~p,t)e l ^-^ q -^ t ] + c.c. 
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where the sum in ^ has been replaced by an integral, since we are interested in the 
thermodynamic limit N = oo. For small q- values, which are, by the way, responsible 
for the asymptotic behavior, the above equation reduces to 

T(q, t) « q 2 dp \C( P , t)C(q - p, f) e 2 W (8) 

J — TT 

+ C{p,t)C*(q-p,t)e-^' (p *> qt ] + c.c. 

where lo'{p) = The first term in the r.h.s. and its complex conjugate are negligible 
because of the rapidly varying phase. By further neglecting the dependence of the 
sound velocity on q (i.e. w'(p) ~ 1), one finds that T(q, t) has the following structure, 

T{q,t) = M(q,t)e lqt + c.c. (9) 

where 

M(q,t) = q 2 ( dpC*(p,t)C(q-p,t) (10) 



is the kernel in the new description. Upon substituting equation ([9]) and ([6]) into ([2]), 
one obtains, in the slowly varying envelope approximation, C <C qC, 

2C(q,t)+e [ dt' M(q,t-t')C(q,t') = (11) 
Jo 

and a similar expression for C* . 

Equation (|TT|) have been obtained after discarding the second order time 
derivative of C(q, t) as well as the integral term proportional to C, besides all rapidly 
rotating terms. The validity of this approximation is related to the separation between 
the decay rate of C(q,t) and u>(q); its correctness will be cheked a posteriori, after 
discussing the scaling behaviour of C(q,t). Notice also that in this approximation, 
Umklapp processes do not contribute: it is in fact well known that they are negligible 
for long-wavelength phonons in ID [TT] . 

Having transformed the second order differential equation into a first order one 
for C(q, t), it is now possible to develop a scaling analysis for the dependence of C on 
q and t, owing to the homogeneous structure of the resulting equation (see also [28] . 
where a similar equation was investigated) . More precisely, one can easily verify that 
the following Ansatz holds, 

C(q,t) = g(^tq 3 / 2 ) , M(q, t) = q a f( ^etq 3 ' 2 ). (12) 

This shows that the decay rate for the evolution of C(q,t) is given by q 3 ^ 2 ^/e. Since 
the corresponding rate for the phase factor is q, one can conclude that amplitude and 
phase dynamics are increasingly separated for q — > 0. High q- values (q w 1) are those 
for which the slowly varying envelope approximation is less accurate. However, if e is 
small enough, such modes are correctly described, too. This has been checked in the 
numerical solution of equation @ (see the following subsection) . 

The functions / and g can be determined by substituting the previous expressions 
into (|TTj) and fTU]) . Upon setting x = y/etq 3 / 2 , one obtains the equation 

Jr = - f dyf(x-y)g(y) (13) 
dx Jo 

r+oo 

f(x)=x-y s dyg*(\x 2 / 3 -y\ 3 / 2 )g(yV 2 ) , (14) 
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where the integral in (fTU|) has be extended to infinity. 

The asymptotic behaviour for x — > can be determined analytically. Actually, 
one can neglect the dependence on x in the integral appearing in (|14|) , thus obtaining 

/(*) = -h ( 15 ) 



where 



which implies that 



+00 

dy\g\ 2 (y 3/2 ) (16) 



2 

M(q,t) = ^^ (17) 

By replacing equation fj 1 5[) into (|13j) and introducing the dummy variable z = x — y, 
one finds 

dg f x , g{x-z) 1/3 3 dg 4/3 

^=-«y o dz- p[ j S - = -3agz/ - -a-x +... (18) 

In the small- a; limit, the last term in the r.h.s. is negligible with respect to the term 
in the l.h.s. After solving the resulting differential equation one obtains, 

g(x) = -exp( —J (19) 

where the multiplicative factor has been determined by imposing the normalization 
condition G(0) = 1. 

Under the approximation that expression (| 1 5() holds for every x, one can exactly 
solve equation (|13|) by Laplace transformation. As a result, one obtains 

iz 1 / 3 

C ^ z ) = ■ 4/3 T 2 • ( 2 °) 

This expression is precisely the Laplace transform of the Mittag-LefHer function 
E^-iXty) [301 ES] for n = 4/3 and A = (a 9 2 ) 3 / 4 . The function £J M interpolates 
between stretched exponential and power-law decays (see Appendix B in |29j). The 
asymptotic form for small values of the function argument yields the same dependence 
for C (see the first of equations ([T2])). This observation suggests that the effective 
evolution of fluctuations should be represented by the fractional differential equation 

§^C(q,t) + A"C(g,i) - 0. (21) 

The case of interest here (1 < fj, < 2) corresponds to the so-called fractional oscillations 
[30] . It should be emphasized that in the present context, memory arises as a genuine 
many-body effect and is not postulated a priori. 



3.2. Numerical analysis 

To assess the validity of the above calculation, we have numerically integrated 
equations ([2|) by the Euler method for the original dispersion relation u)(q) and different 
values of e. We have verified that a time step At = 0.01 guarantees a good numerical 
accuracy over the explored time range. The Fourier transform G(q, u>) is plotted in 
figure [1] for three different q- values versus lo — wjvf(g), where u>m{i) is the frequency 
corresponding to the maximal value Gm of the spectrum (this is equivalent to removing 
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Figure 1. Fourier transform G(q,u>) of the correlation functions for e = 1, 
N = 2000 and the three wavenumbers, q = 7r/100 (circles), 2-7r/100 (squares), 
and 37r/100 (diamonds). The solid line corresponds to the approximate analytic 
theory, i.e. it is obtained by Fourier transforming the function C(q,t) defined in 
I|12ll9j l. The same curves are plotted in the inset in log-log scales, where only 
positive frequencies are shown. 



the oscillating component from G(q,t)). Furthermore, in order to test relation Q12]). 
the vertical axis is scaled to the maximum G-value, while the frequencies are divided 
by the half-width j(q) at half of the maximum height. This latter quantity can 
be interpreted as the inverse lifetime of fluctuations of wavenumber q. The good 
data collapse confirms the existence of a scaling regime. Moreover, the numerical 
data is compared with the approximate analytical expression, obtained by Fourier 
transforming the function defined in (fT2")l and (119[) (solid curve) . The overall agreement 
is excellent; there are only minor deviations for small values of to — com where (|19p 
is not strictly applicable. Moreover, in the inset of figure [TJ where the same curves 
are plotted using doubly logarithmic scales, one can see that the lineshapes follow the 
predicted power-law, a;" 7 / 3 , over a wide range of frequencies. This scaling has been 
deduced by Fourier transforming g(x) in equation (|19[) . 

In figure [2] we show that 7(g) is indeed proportional to sfeq 3 ^ 2 . It is noteworthy 
to remark that the agreement is very good also also for a relatively large value of the 
coupling constant such as e = I , although the slowly varying envelope approximation is 
not correct for large wavenumbers. The deviations observed at small q- values for small 
couplings are due to the very slow convergence in time. Better performances could be 
obtained by increasing both N and the integration time (10 , in our units), however, 
well beyond our current capabilities. Finally, we have plotted the Fourier tranform 
of the memory function (see figure [3]) to check its scaling behaviour. By Fourier 
transforming equation (| 1T|) . one finds a behaviour like uo^ 1 / 3 which is reasonably well 
reproduced in figure [3]). Once again, in order to obtain a larger scaling range, it would 
be necessary to consider significantly larger system sizes. 

Making use of equations (|12ll9p . it can be shown that the memory function T 
contains terms of the form q 2 e ±lqt jt 2 ! 3 , i.e. it oscillates with a power-law envelope. 
Accordingly, its Laplace trasform has branch-cut singularities of the form q 2 /(z±q)3 . 
This finding is not consistent with the heuristic assumption made in [T4J, [15] and 
the result of [17], where simplified MCT equations were solved with the Ansatz 




Figure 2. Scaling of the linewidth 7(g) of G(q, 10) with q for three different values 
of the coupling constant e and N = 2000. The solid line corresponding to the 
expected power law q 3 ' 2 is plotted for reference. 




Figure 3. Fourier transform of the memory function T(q,w) for the same 
parameter values as in figure [l] and for three wavenumbers, namely, 7r/100 (solid 
line), 2-71-/IOO (dashed line), and 37r/100 (dotted line). The ui axis is suitably 
shifted in order to get rid of the oscillatory peak at frequency q, so that T(q,w) 
reduces to M(q,u>). For comparison, the expected power-law 

i s 

also 

presented (thin solid line). 



T(q, z) — q 2 V{z). In addition, the numerical solution does not show any signature 
of the q 2 /z 1 / 3 dependence found in [17]. For instance, it would imply a peak at lj = 
in the spectrum of T which is, instead, absent in our data. 

4. Quartic nonlinearities 

In this section, we turn our attention to the case of a vanishing 53. This the case of 
symmetric potentials (with respect to the equilibrium position) , but it is more general 
than that. The kernel of the mode coupling equations writes in this case 

T(q,t)= u 2 (q)(^) 2 J2 G(p,t)G(pi,t)G( P ",t) . (22) 

p+p' +p" — q— 0, ±7r 



Anomalous kinetics from ID MCT 9 

Equation ^ is still formally valid as it stands, with the (bare) coupling constant now 
given by e = 15(fl 4 k B T /2tt) 2 [15]. 

Analytical results 

One can repeat the same analysis performed for the cubic case. By substituting 
expression ([6]) into (|22[> . we must now deal with a double integral, 

|>7T |>7T 

T(q,t) =tu 2 {q) dp dp' C{p,t)C{p\t)C{q-p-p\t)e l ^ p '> +UJ ( p '') + ^'i- p - p 'K t 



+C{p, t)C{p', t)C*(q -p-p', i )e i ^(p)+"(f')-^(9-p-p')]i 
+C(p, t)C*(p', t)C(q -p-p', f) e *Kp)-"(p')+"(9-P-p')]f 
+C(p, t)C* 0', t)C*(q-p- P ', t yl"{p)-"{p')-"(i-p-p')]t + c.c 1 (23) 

By neglecting the non-resonant terms in the limit g — > 0, the above expression 
simplifies to 

T(q,t) =q 2 P dp P dp'\c*{p,t)C*{p',t)C{q-p-p\t)e^^ +p '^ + c.c] (24) 



which, by recalling that u' ~ 1, is formally equivalent to equation ©. Altogether, 
equation pip applies also to the quartic potential, the kernel M being now given by 

M{q,t) = q 2 j dp f dp'C*(p,t)C*(p,t)C(q-p-p',t) (25) 

Equations (|11I25|) satisfy the new scaling relations 

C(q,t) = g{V~etq 2 ) , M(q,t) = q 4 f(V^tq 2 ) (26) 

As in the cubic case, the functions / and g can be determined by substituting the 
previous expressions into equations (|11I25|) . Upon setting x = \feq 2 t one obtains, 

^ = -J*dyf(x-y)g(y) (27) 

/+oo r+oo 
dy dy'g*{y 2 )g*{y' 2 )g{\x^-y-y'\ 2 ) (28) 

-oo J — oo 

In the limit x — > 0, one can analytically derive the asymptotic behaviour that is simply 

f{x) = - (29) 
x 

where 

f+oo r+oc 

dy / dy'g*(y 2 )g*(y' 2 )g(\y + y'\ 2 ) (30) 
which, in turn, implies that 

%*) = ^y (31) 

By inserting equation (j2"9"]) into (|2T|) one finds 

-7- = - f dya-^— = -a[ dy giy)S(x - y) = -ag(x) (32) 
dx J x-y J Q 
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whose, properly normalized solution, is 



g{x) = - cxp(-aa:) (33) 



By assuming that (|29[) is valid in the entire x-range, we find that the Laplace transform 
of C(q, t) is 

iz + bq z 

where 6 is a suitable constant. Remarkably, the correlation function is a simple 
exponential, and its Fourier transform is a Lorentzian curve. Note that also in this 
case we obtain the Mittag-Lefner function with index fi = 1, which is precisely an 
exponential [50] , 



4-2. Numerical analysis 

Also for the quartic case we have compared the analytical solutions with the numerical 
integration of the full mode-coupling equations. In figure|3]we have plotted the Fourier 
transform of the correlation function G(q,uj), by following the same approach used in 
the cubic case, for three different g-values. The y axis is scaled by Gm and the 
frequencies arc divided by half the maximum width 7(g). The good data collapse 
confirms the existence of a scaling regime (see (|26p ). Moreover, at variance with the 
previous case, the approximate analytic solution agrees very well with the numerical 
data in the entire range, including the region around the maximum, where correctness 
is not a priori obvious. 

A further test of the theory can be made by determining the behavior of the 
linewidths 7(9), which our theory predicts to scale as \feq 2 ■ In order to verify this 
scaling behaviour, it is necessary to test small wavenumbers, but the smaller q, the 
longer it takes to reach the asymptotic regime. Moreover, the rate of convergence is 
controlled also by the coupling strength e. In figure [SJ we plot the data corresponding 
to four different values of e. They show that already for e = 0.1, the chosen time 
span (10 4 units) is not enough to ensure a reasonable convergence for the smaller 
g-values. Nevertheless, a fit of the tail of the curve corresponding to e — 1 reveals 
a reasonable agreement with the theoretical prediction (the slope is indeed equal to 
1.91). Unfortunately, choosing yet larger values of e would descrease the width the 
scaling range. 

Finally, in figure [5] we plot T(q,uj) versus u> — ujm . At variance with the cubic 
case, the spectrum does not exhibit any low-frequency divergence. This is consistent 
with the analytic solution, since the Fourier transform of equation (|3ip saturates for 
uj -> 0. 



5. Energy current correlator 

In the previous sections we have solved the mode-coupling equations as it applies to 
monoatomic chain of anharmonic oscillators for two distinct leading nonlincarities. 
Here, we use such results to determine the transport properties of the underlying 
lattice. In particular, we are interested in the thermal conductivity as defined by the 
Green-Kubo formula, 

POO 

k oc / (J(t)J(0))dt (35) 
Jo 
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Figure 4. Quartic case: Fourier transform G(q,uj) of the correlation function 
for e = 1, N = 2000 and the three wavenumbers q = 7r/100 (circles), 27r/100 
(squares), and 3-7r/100 (diamonds). The solid line corresponds to the approximate 
analytical solution obtained by Fourier transforming the expression defined in 




Figure 5. Scaling of the linewidth f(q) of G(q,Lo) for four different values of the 
coupling constant e and N = 100. The curves are obtained by integrating over 
10 4 time units. The solid line with slope 1.91 is a fit of the data corresponding to 
e = 1. 



where J{t) is the total heat current, and (. . .) brackets denote the equilibrium 
average. A general expression for the heat current has been derived in [11] , For the 
determination of the scaling behaviour it is sufficient to consider only the harmonic 
part, which can be expressed as a sum over wavevectors 

J = Y,b(q)Q(q,t)P*(q,t) (36) 
i 

where Q(q,t), and P(q,t) — Q{q,t) are the canonical variables in Fourier space and 
b{q) = iu>{q)^- . (37) 

This amounts to disregarding higher-order terms which are believed not to modify the 
leading behaviour. Under the same approximations that allow deriving the mode- 
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Figure 6. Fourier transforms of the memory function for three values of the 
wave-vector for q = 7r/100 (solid line), 27r/100 (dashed line), and 37r/100 (dotted 
line). For small frequencies they tend to a constant as we expect by analytic 
solution. 



coupling equations, i.e. by neglecting correlations of order higher than two, one 
obtains, 

(J(t)J(0)) = I %) I 2 {(Q(q,t)Q(q,0))(P(q,t)P*(q,0)) + (Q(q 7 t)P* (q,0)) 2 } (38) 

This expression can be further simplified under the assumption P(q) ~ uj(q)Q(q), that 
is certainly valid in the small g-limit. Altogether, this leads to the expression proposed 
in [17], 

(J(t)J(O)) « J2 (t)^ 2 ^) (39) 
<; 

By taking into account the general expression ([6]), one obtains 

(J(i)J(O)) cx I ' dq\c 2 { qi t)e 2tqt +C* 2 (q,t)e- 2lqt +2C(q,t)C*(q,t)] (40) 



where, as usual in the thermodynamic limit, we have replaced the sum with an integral 
and have made use of the fact that for small q (the integral is dominated by the small-q 
terms), uj(q) « q. We can finally use (|12p to determine the scaling behaviour of the 
heat flux autocorrelation. Since the contribution of the terms proportional to e 2 * 9 * is 
negligible, due its "rapid" oscillations, one obtains, 

(J(i)J(O)) CX J dqg 2 (V^tq 3 / 2 ) cx ^ (41) 

The resulting integral in equation (J3SJ) thereby diverges, revealing that the effective 
conductivity of an infinite lattice is infinite. In order to estimate the dependence of k 
on the system size L it is customary to restrict the integration range to times smaller 
than L/c where c is the sound velocity. As a result, one finds that k(L) diverges as 
L 1 / 3 , i.e. a — 1/3. Therefore we conclude that mode coupling predictions are now in 
full agreement with those of the renormalization group [12] . 

A similar analysis for the quartic case, leads, using equation ([26]) . to 



(J(t)J(O)) cx J dqg 2 (q 2 t)K^ (42) 
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Figure 7. Correlation function of the heat-flux for the cubic (a) and quartic (b) 
case. The curves have been obtained for TV = 2000. In (a), also the outcome of 
a calculation made with N = 4000 are plotted to show the good overlap in this 
time interval; the major difference is the obvious shift of the peak observed for a 
time equal to the system length divided by the sound velocity. 



which, in turn, corrsponds to k ~ L 1 / 2 , i.e. a = 1/2. Finally, we have compared 
the theoretical predictions with the computation of the sum (|39|) complemented by 
direct integration of the mode coupling equations for G(q, t) (see figure[7]). A fit of the 
scaling region yields a slope equal to —0.69 in (a) and —0.55 in (b) to be compared 
with the theoretical expections that are respectively equal to —2/3 and —1/2. The 
deviations are analogous to those ones observed in the previous explicit analysis of the 
behaviour of G(q,t) itself and can thus be attributed to finite size/time corrections. 

6. Conclusions 

In this paper we solved exactly the ID MCT corresponding to a chain of nonlinear 
oscillators. The only assumption we made is that the dynamics of both the correlator 
G(q,t) and the memory function T(q,t) factorizes in the product of a "rapidly" 
oscillating term times a slowly decaying function, see equation ([6]). This Ansatz 
is confirmed a posteriori by the scaling of the analytic solution and by numerical 
calculations. The oscillating behaviour of T(q, t) represents an element of novelty that 
was disregarded in previous studies [TSl E] (see hovever the perturbative calculation 
in [H]). In fact, taking the oscillations into account, we find that the decay of the 
envelope M(q,t) depends on the order of the leading nonlinearity of the theory. In 
the cubic case, M oc £~ 2 / 3 , while a faster decay, M oc f -1 , is found in the quartic 
one, see equations (|17I31|) . As a consequence, in the two cases, the decay of G(q, t) is, 
respectively, anomalous (i.e. nonexponential) and exponential, see equations (|19I33|) . 
Upon comparing with the solution of the Markovian equation ([3]) , we thus argue that 
the cubic case is characterized by a diverging viscosity 77, while a normal hydrodynamic 
behaviour (finite rf) is found in the quartic case. At the same time, the other transport 
coefficient, the thermal conductivity «, diverges in both cases. Indeed, the energy- 
current autocorrelation (as estimated from equation ([39]) ) shows a long-time tail, t~ 2 / 3 
and t~ x / 2 in the cubic an the quartic case, respectively. This implies that for large L, 
k diverges as L Q , with a = 1/3 (cubic nonlinearity) or a = 1/2 (quartic nonlinearity). 
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Since generically, the leading nonlinear term is cubic, we can conclude that our 
analysis reconciles MCT with the renormalization-group prediction PJ, ■ Accordingly, 
on the theoretical side there is no longer any doubt on the existence of a wide 
universality class characterized by the exponent a = 1/3. 

The quartic case is less generic, since it requires the fine tuning of a control 
parameter (as for phase transitions) or the presence of a symmetry forbidding the 
existence of e.g. odd terms in the potential. This is, for instance, the case of the 
FPU-/3 model which, by construction, contains only quadratic and quartic terms. In 
fact, the most accurate numerical findings obtained for the FPU-/3 model [5T] are 
qualitatively compatible with the scenario emerging from the analytic solution of the 
MCT in the quartic case, since a > 2/5 and does not exhibit a clear saturation with 
the system size. 

The two-universality-class scenario emerging from our solution of the MCT is 
consistent with that one outlined in [3T]. There, the scaling behavior has been 
estimated by complementing hydrodynamic and thermodynamic arguments under 
suitable assumptions. Such an analysis, that applies to quartic nonlinearity, predicts 
the same behavior for k as ours, as well as a finite rj. As the Authors of reference 
[31] pointed out, the cubic case requires a fully self-consistent calculation which is 
precisely the task we accomplished here. Moreover, as already remarked, we conclude 
that both transport coefficients diverge. In addition, we argue that since M and 
(J(i)J(O)) decay with the same power-law, 77 should scale in the same way as K. A 
few words are also needed to comment about the definition of the second class, which 
is identified in the present paper from the vanishing of the cubic term in the expansion 
of the effective potential {T]), while in [31] from the equality between constant-pressure 
and constant-volume specific heats. At least in the FPU-/3 model both constraints 
are simultaneously fulfilled. More in general, it seems plausible that these are two 
equivalent definitions of the same class but it is not certain. 

Altogether, one of the important consequences of our analytical treatment is that 
fluctuations and response of ID models cannot be accounted for by simple equations 
such as ([3]). Indeed, memory effects emerging from the nonlinear interaction of 
long-wavelength modes are essential and cannot be considered as small corrections. 
As a consequence, ([3]) should be replaced by equations like ([2T]) . meaning that 
hydrodynamic fluctuations must be described by generalized Langevin equation with 
power-law memory [32j [33] . This provides a sound basis for the connection between 
anomalous transport and kinetics (superdiffusion in this case) [34] . 




Finally, a few words on the consistency between this theoretical analysis and 
the simulation of specific models. Over the years it has become clear that numerical 
studies are affected by various types of finite-size corrections which make it difficult if 
not impossible to draw convincing conclusions (and this is one of the main reasons why 
we turned our attention to analytical arguments). The ubiquity of slow corrections 
suggests to carefully consider all the implicit assumptions made in the various theories, 
in the hope to understand the possibly universal origin of such deviation or even 
(though unlikely) to obtain corrections to the leading exponents. In this perspective, 
we wish to point out that the MCT solved in this paper is based on the assumption that 
the hydrodynamic properties are dominated by the coupling among sound-modes. A 
more general and complete mode-coupling analysis should consider also the coupling 
with thermal modes. General considerations indicate that the asymptotic scaling 
properties are not affected by lower-order corrections emerging in this more general 
approach. On the other hand, numerical experiments seem to suggest that such 
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corrections could affect significantly finite size/time corrections to the scaling. 
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